Dynamic Scaling of Bred Vectors in Chaotic Extended Systems 
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We argue that the spatiotemporal dynamics of bred vectors in chaotic extended systems are 
related to a kinetic roughening process in the Kardar-Parisi-Zhang universality class. This implies 
that there exists a characteristic length scale corresponding to the typical extend over which the 
finite-size perturbation is actually correlated in space. This can be used as a quantitative parameter 
to characterize the degree of projection of the bred vectors into the dynamical attractor. 
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A standard tool for studying chaotic behavior in dy- 
namical systems is the computation of the characteristic 
Lyapunov exponents, which, roughly speaking, measure 
the typical exponential growth rate of an infinitesimal 
disturbance |2| ■ The characteristic Lyapunov expo- 
nents in extended systems are defined in a similar way 
as their low-dimensional counterpart and can be calcu- 
lated from the linearization of the equations of motion 
I^Q. The main point is that the growth of an infinites- 
imal perturbation is described by the linear equations 
for the tangent space, the so-called Lyapunov vectors. 
However, for many practical purposes, Lyapunov vec- 
tors may be irrelevant as indicators of, for instance, the 
predictability time. Indeed, in realistic situations, the 
error in the initial condition is finite. The important 
fact is that the evolution of finite errors is not con- 
fined to the tangent space, as defined by the growth of 
linearized perturbations, but is controlled by the com- 
plete non-linear dynamics. In order to deal with real- 
istic perturbations, the concept of finite size Lyapunov 
exponents has been found 5] to be useful to analyze 
predictability in high-dimensional systems 3] . 

A good example with important practical application 
occurs in weather forecasting research, where one deals 
with the whole Earth's atmosphere- an extremely high 
dimensional system in which initial conditions can be 
determined only with limited accuracy. In this context, 
Patil et. al. have recently introduced the concept of 
hred vector (BV) jfl, (Tj , as the spatio-temporal evolu- 
tion of a statistical ensemble of non-infinitesimal per- 
turbations. From the analysis of real data provided by 
the US National Weather Service, Patil et. al. have 
proposed and measured an effective local finite-time di- 
mension, by means of the so-called local BV dimension 
statistic 0, . Regions of low BV dimension are iden- 
tified as more predictable than locations where dimen- 
sion is high. The method has recently been applied to 
2D coupled map lattices in which nonlinear time se- 
ries analysis has been used to make local predictions of 
trajectories at low dimension sites 0- Breeding tech- 
niques constitute an important tool in modern weather 
forecasting research and one can envisage that BVs 
can be used in other spatially extended systems with 



chaotic dynamics where local short time forecasts can 
also be feasible j^. However, many questions concern- 
ing statistical and dynamical properties of BVs are still 
unanswered. The aim of this Letter is to shed some 
light into these questions. 

In this Letter we study the spatiotemporal dynam- 
ics of homogeneous finite-size errors and focuss on the 
propagation dynamics of BVs in chaotic extended sys- 
tems. We argue that, after a suitable transformation 
of variables, the BV dynamics can be interpreted as a 
kinetic roughening^rocess in the Kardar-Parisi-Zhang 
universality class |9(. We also find that the breeding 
procedure introduces a characteristic length scale cor- 
responding to the typical extend over which the finite- 
size perturbation is actually correlated in space. This 
can be used as a quantitative parameter to characterize 
the degree of projection of the BVs into the dynami- 
cal attractor, which is of major interest in probabilistic 
forecasting techniques based on ensembles of BVs. 

We exemplify our results by means of numerical simu- 
lations of coupled map lattices in one dimension, which 
are simple model systems exhibiting space-time chaos 
and convenient as far as the computing time is con- 
cerned. Then, we shall be considering a coupled map 
array consisting of L chaotic oscillators given by 



u{x, t + 1) 
vf{u{x-l,t)) - 



uf{u{x + l,t)) + 
{l-2v)f{u{x,t)) 



(1) 



where x = 1, 2, . . . , L, f{u) is a chaotic map, u is the 
coupling constant, and periodic boundary conditions 
are imposed. We have fixed the coupling to v = 1/3 
in all the simulations presented in this Letter. We 
have carried out simulations for two different choices 
of the map, the chaotic logistic map f{u) = 4u(l — u), 
< u < 1 and the tent map f{u) = 1 — 2\u — 1/2|, 
< u < 1. For the sake of brevity, all the results we 
present below correspond to coupled logistic maps, but 
similar results were obtained for the tent map. 

Scaling of finite perturbations.- Let us first consider 
the evolution of finite size homogeneously perturbed 
trajectories in our model system (Q. Given an ini- 
tial condition u°(a;,0), the solution u'^{x,t) is univo- 
cally determined by computing Eq.l^ for a number 
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FIG. 1: Numerical results for the propagation of finite-size 
errors in the coupled logistic maps. Upper panel shows the 
amplitude factor e{t) vs. time for perturbations starting 
with initial amplitudes of eo = 10"^, 10"^ 10"'', 10"^, 10"" 
and 10"^^ (from top to bottom) in ID lattices of L = 1024 
sites and results were averaged over 600 different initial con- 
ditions. Lower panel shows the variance cr^(t) for the same 
initial perturbations as before and eo decreasing from top 
to bottom. The same data are plotted in log-log scale in 
the inset to show the power-law behavior (j'^{t) ~ t'^'^ . The 
straight line has a slope 0.60 and is plotted to guide the eye. 



t of time steps. This will be our reference trajectory 
and we shall be studying the evolution of finite per- 
turbations around that reference solution. Since we 
are interested here in the propagation of real (non- 
infinitesimal) errors, we should avoid linearization of 
Eq.Q. Instead, we compute the trajectories gener- 
ated by iterating for an ensemble of initial condi- 
tions u{x, 0) — 0) -f Su{x, 0), where 6u{x, 0) is uni- 
formly distributed in (— eo,eo). For each iteration of 
the lattice ^ the difference 5u(x, t) = u{x, t) — u'^{x, t) 
between the reference trajectory and every one of the 
disturbed solutions is calculated and represents the evo- 
lution of homogeneously distributed finite errors in the 
system. Although the disturbances are initially homo- 
geneous and uncorrelated in space, as time goes by, they 
propagate, get correlated, and grow exponentially in 
size. The statistical fluctuations of disturbances can be 
characterized by studying now the ensemble of finite 
perturbations {(5u„(a:, t)}^^j^, which correspond to N 
independent realizations of the initial perturbation. In 
our simple model we can use ensembles of hundred of 
samples to obtain a good statistics. 

Firstly, we introduce what we call the amplitude fac- 
tor e(t) as the spatial geometrical mean value of the 
perturbation: 



e{t)^ll\Su{x,t)\ 



l/L 



(2) 



which, as we shall see below, contains the information 
about the dominant exponential growth rate. In Fig. 1 
we plot loge(i) vs. time for different values of the initial 
perturbation amplitude cq. We can see that for times 
t < iiin(eo) the average amplitude factor grows exponen- 



tially in time €{t) « eoexp(At). We demonstrate below 
that A indeed corresponds to the maximal Lyapunov ex- 
ponent. For longer times, t > tiin(eo) the amplitude fac- 
tor saturates to a constant value. Both, the saturation 
constant and the maximal Lyapunov exponent, are in- 
dependent of the initial perturbation size eg. However, 
the saturation times iiin(eo) increase as the size of the 
initial perturbation eo becomes smaller. This indicates 
that iiin corresponds to the crossover time at which the 
dynamics of a finite size perturbation depart from the 
linear approximation {i.e. the Lyapunov vectors) . This 
crossover occurs because the tangent Lyapunov vectors 
describe only the behavior of strictly infinitesimal per- 
turbations (eo 0). Therefore, for times t > tiin(eo) 
nonlinear corrections, due to finiteness of the initial per- 
turbation, come into play and force the errors out of the 
tangent space. From then on, the linear approximation 
cannot describe the evolution of errors. One then ex- 
pects that tiin(eo) ^ oo as eo ^ 0. 

Besides exponential growth, correlations are dynam- 
ically generated during the evolution of perturbations. 
Correlations contain information about the sub-leading 
Lyapunov exponents and thus also contribute to the 
perturbation size growth. The important role of corre- 
lations can be better realized after subtraction of the 
dominant exponential growth component given by e{t). 
We find that a very useful indicator is given by the re- 
duced perturbations Sr{x, t) that we define as 



6r{x, t) 



5u{x, t) 



(3) 



where the dominant exponential growth is globally re- 
moved and one is left with the effect of correlations. 
Statistical fiuctuations of the reduced perturbations are 
measured by the variance cr^(t) = (i5r(a;, i)^), where 
(• • • ) stands for average over realizations of the initial 
perturbation and the over bar is a spatial average. As 
we shall see below, a{t) gives information about the 
growth of perturbations due solely to correlations. 

In Fig. 1 we show our numerical results for the vari- 
ance of the reduced perturbations, Eq.©. The vari- 
ance (T^(i) grows as an exponential power law (T^(i) ~ 
exp{t'^^) for times t < tiin, i.e. before saturation by 
finiteness of the initial disturbance. A least-squared 
fit of log(log((T)) vs. log(i) gives the exponent /3 = 
0.30 ± 0.05 (see the inset of Fig. 1). As before, this 
rapid growth occurs for times such that the dynamics 
of disturbances are well described by the Lyapunov vec- 
tors. 

The magnitude and extent of spatial correlations 
can be measured by use of the site-site correlation 
function p{x,t) — {Su{xo,t)Su{x -\- xo,t))/{Su{xo,t)'^). 
The magnitude of spatial correlations increases in time 
p(x,ti) < p{x,t2) if ti < t2 for times ti,t2 < tun. How- 
ever, correlations become progressively smaller for times 
t > tiin . This indicates that the effect of having initially 
finite perturbations is to introduce a characteristic time 
Uini^o) marking the typical time it takes for the system 
to depart from tangent space (with the building-up of 
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correlations) to truly non-linear evolution (uncorrelated 
errors). For times larger than tun perturbations quickly 
get spatially uncorrelated, as can be seen from the de- 
cay of the reduced variance in Fig. 1. Further analysis 
[lOj shows that errors become undistinguishable from 
actual white noise. 

Our numerical results can be explained analytically 
by making use of a transformation first proposed by 
Pikovsky and Politi for the actual Lyapunov vectors 
[ill fl^ . We find that the dynamics of finite perturba- 
tions can also be seen as a kinetic roughening processes 
of the surface defined by h{x,t) = \og\Su{x,t)\. This 
is a very useful transformation that allows us to make 
use of existing results in the field of nonequilibrium sur- 
face growth. Indeed, Politi and Pikovsky 0, have 
shown that errors in many extended systems lead to 
surface growth process in the universality class of KPZ 
0. Let us first consider the amplitude factor defined 
in Eq. ij^Jl and show that it is related to the average 
surface velocity. We can write \5u{x,t)\ = exp[h{x,t)] 
and thus e{t) = exp[(l/L) /i(a:, t)] = exp[/i(t)]. 

We then obtain that the amplitude factor must grow as 
e{t) = Co exp(At), where A is the largest Lyapunov expo- 
nent and corresponds to the surface velocity, in agree- 
ment with our numerical results shown in Fig. 1. 

Also the time behavior of the variance of the re- 
duced perturbations shown in Fig. 1 can be analyti- 
cally related to surface scaling properties. We find that 
(T(t)2 = V2^exp[2W{L,t)^], where W{L,t) is the sur- 
face width and depends on the system size L. This 
result is easily obtained by assuming a Gaussian distri- 
bution of surface heights. KPZ behavior then implies 
that the width scales a,s W{L,t) ~ for times t < ts{L) 
and saturates to a size dependent value, W{L, t) ~ 
for t > ts{L), where /3 and a are the growth and rough- 
ness exponent respectively that are known to have the 
values (3 = 1/3 and a = 1/2 for the KPZ universality 
class in one dimension 0. This explains the exp{t'^^^) 
scaling observed in the inset of Fig. 1. It is worth men- 
tioning that all the numerical results presented here are 
obtained for system sizes such that inn ^ tg, so that 
the existence of could be clearly seen. 

Dynamic scaling of bred vectors.- An important 
problem where finite initial errors are a concern is in 
weather forecasting. In this context, initial configura- 
tions that belong or are very close to the dynamical 
attractor are much sought after. These prepared config- 
urations are then used as initial conditions in the atmo- 
spheric models in order to have a probabilistic forecast 
for a certain time window. Predictions are expected 
to be better when the initial conditions are closer to 
the dynamical attractor of the system. To achieve the 
goal of a maximal projection, the computer runs of the 
atmospheric models are started from random pertur- 
bations and BV techniques are used to allow the 
perturbed trajectories to get closer to the chaotic at- 
tractor. BVs are defined in analogy to the operation of 
data reassimilation in numerical models of atmospheric 
evolution, in which the output of the numerical model is 



corrected by the observed experimental data after short 
periods of time {ti, T2, T3..}. BVs are defined by multi- 
plying each member of the ensemble by a reduction fac- 
tor {fci, A;2, fcs---} at those times so that the amplitude 
factor of the perturbed solutions eo is kept constant. 
Several criteria could be chosen to construct BVs. We 
define here what we call continuous bred vectors 5B{x, t) 
as the perturbation obtained by dividing by the am- 
plitude factor at every time step. Since the dominant 
exponential growth is filtered out, BVs evolve with a 
constant amplitude factor by definition, but exhibit the 
actual spatial correlations (see for further details). 

We have studied the dynamics of BVs in the coupled- 
map lattice model given by Eq. As in the case 
of finite perturbations we have considered a reference 
trajectory vP{x,t) and its corresponding ensemble of 
randomly perturbed trajectories u{x,0) = u*'(a;,0) 4- 
(5it(a;,0). But now breeding is applied so that we define 
the re-scaled error (the BV) 5B{x,t) = Su{x,t)/e{t) at 
every time step and the perturbed trajectories are then 
evolved by introducing u'^{x,t) + €o5B{x,t) in Eq. ^ 
to obtain u{x,t + 1). We have carried out simulations 
of the model for different values of the finite amplitude 
eo, which is the only external parameter. The imple- 
mented procedure is intended to mimic the one used in 
probabilistic weather forecasting 0. 

The reduced perturbations 6r{x, t) = eo dB{x, t) and, 
in particular its variance, give a direct measure of the 
spatial fluctuations of the BVs. As we have shown 
in our previous analysis of finite perturbations, which 
also applies here, the variance cr^ {t) grows with time as 
<7^ ^ exp{t'^^) until it saturates due to nonlinear terms 
becoming important in the equation of evolution of fi- 
nite errors. The main point now is that, in contrast with 
the case of bare finite-size perturbations, the breeding 
procedure leads to a non-trivial stationary value of the 
statistical fluctuations of BV. More importantly, this 
stationary state is characterized by having space corre- 
lations, as can be seen from the finite stationary values 
of the variance tr^ in Fig. 2. The time scale r(eo,i) 
characterizes the crossover to the stationary regime, 
where fluctuations are no longer time dependent. After 
that time BVs are correlated over regions of character- 
istic size l{eo, L) ^ t^^^ , where z = a/(3 is the dynamic 
exponent. We find that numerical data of BV fluctua- 
tions can be cast in a dynamic scaling ansatz given by 

\og[<j\t,eo)]=t'f'g{t/T{eo)), (4) 

where the scaling function g{v) ~ v~'^^ ii v ^ 1, and 
G{v) ~ const if w <C 1. The characteristic time is found 
to scale as T(eo) ~ | log(eo)|^/'', as shown by collapse of 
the data in the inset of Fig. 2. 

We can now make use of the surface growth map- 
ping, h{x,t) = log\SB{x,t)\, to obtain some analyti- 
cal understanding on the dynamic scaling behavior of 
BVs. Among other quantities, one can obtain the scal- 
ing law for the saturation time T(eo,i) as follows. Let 
AI be the upper bound for the linear approximation to 



4 




logio(t) 

FIG. 2; Dynamic scaling of BVs for the coupled lo- 
gistic maps. Main panel shows the variance u^{t) vs. 
time for perturbations kept at finite amplitudes eo — 
10"^10"'',10"^10"^ and 10"^ (from bottom to top) by 
means of the breeding procedure as described in the text. 
The straight line has a slope 0.60 and is a guide to the eye. 
The inset shows a data collapse in log-log scale according 
to Eq. @ with exponent f3 = 0.28 ± 0.05. Results were ob- 
tained in systems of linear size L = 1024 and averages over 
500 different initial conditions were taken. 

be valid to describe the evolution of perturbations in 
the system. To be precise, if the error is larger than 
M, \5u\ > M, then the time evolution of 5u{x,t) is 
not longer in the tangent space (infinitesimal), but fi- 
nite and to be described by higher order nonlinearities. 
Any finite-size perturbation in the linear regime is then 
bounded by eo SB < M. This implies an upper bound 
in the value that the surface height can take h„iax ^ 
|log(Af/eo)| and therefore a saturation time r corre- 
sponding with the typical time for which the surface 
width W{t) ^ hmax- This leads to a characteristic time 
scale T(eo) | log(Af/eo)|^/'^. This time scale together 
with the usual saturation time Tc ~ associated with 
the kinetic roughening of the surface h are the two char- 
acteristic time scales in the system. Saturation of the 
surface fluctuations is then controlled by the shortest of 
the two, r(eo,i) min{| log(M/eo)p/'3, L^}. This the- 
oretical argument agrees very well with our numerical 
data as can be seen in the inset of Fig. 2, where an ex- 
cellent data collapse is obtained by using our theoretical 
expression for the saturation time r '-^ \ log(Af/eo)|^/'^. 
For clarity of presentation all the values of eo shown in 
Fig. 2 satisfy | log(M/eo)|^/^ < so that saturation 
time is governed by | log(Af/eo)|^^'^. On the contrary if 
smaller systems are used one can see that T \m. 

Due to the scale-invariant dynamics of the sur- 
face h{x,t) = log \SB{x,t)\ there exists a length scale 
l{eo,L) ~ T{eo,Ly^^ that corresponds to the typical 
extend of spatial correlations in the stationary regime. 
This length scale is a quantity of great importance since 



it can be used as an indicator of the degree of projec- 
tion of the BVs into the dynamical attractor. If / < L 
correlations are short ranged and the perturbed trajec- 
tory is only partially projected over the attractor. On 
the contrary, ior I > L BVs are spatially correlated over 
whole system and projection is maximal. In this case, 
BVs dynamics are confined to the tangent space. 

Finally, we can draw some interesting conclusions for 
the generation and application of BVs as ensembles 
of initial conditions in probabilistic forecasting. For a 
given system size L there is a threshold eth ~ exp(L") 
of the fixed amplitude of BVs such that disturbances 
with amplitudes cq < £th propagate as infinitesimal dis- 
turbances, which dynamics are fully described by the 
tangent space equations (Lyapunov vectors). This in 
turn means that BVs are then maximally projected and 
correlated over whole system. In contrast, BVs with 
eo > ^th are non infinitesimal disturbances that cannot 
be described by tangent space equations, but by the full 
nonlinear dynamics. In this case projection into the at- 
tractor is only partial and correlations have a charac- 
teristic size l{eo,L). Our results show that when using 
breeding techniques in probabilistic forecasting for the 
preparation of ensembles of initial conditions, the elec- 
tion of the normalization condition of the BVs is of vital 
importance since it introduces a previously unforeseen 
time/length scale related to degree of projection into 
the dynamical attractor. 
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